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I.  INTRODUCTION 


A.  BACKGROUND 

The  prediction  of  a  ship's  maneuverability  as  a  part  of  the  design  process  is  becoming 
a  reality.  Valuable  feedback  can  be  given  to  a  design  team  during  the  design  process  to 
provide  insight  regarding  a  proposed  design's  performance.  Previously,  the  primary  method 
of  testing  a  hull  design  involved  a  scale  model  subjected  to  various  tests  in  a  wave  tank. 
Construction  and  testing  of  models  can  be  expensive  and  involve  many  man  hours.  It  is 
important  that  once  a  hull's  design  reaches  this  phase  that  the  design  be  as  complete  as 
possible.  If  a  computer  model  is  available  for  simulating  a  proposed  hull's  performance,  the 
design  team's  productivity  can  be  increased.  The  result  will  be  a  more  accurate  scale  model 
for  tank  testing  and,  ultimately,  a  better  ship  design  at  lower  cost. 

Toward  this  end,  many  studies  have  been  conducted  in  various  formats,  with  varying 
degrees  of  complexity.  A  low  order,  easy  to  reproduce,  simplified  model,  coupled  with  a 
scenario  to  test  a  hull's  design,  would  prove  useful  to  design  teams.  This  study  is  aimed  at 
producing  a  simplified  model  of  ship  and  power  plant  dynamics  to  be  used  in  simulation 
studies  of  ship  maneuverability  and  station-keeping.  This  simulation  model  will  then  be  used 
to  study  a  station-keeping  control  algorithm  in  the  presence  of  adverse  disturbances  due  to 
wind  and  wave  conditions.  Further,  this  study  is  aimed  at  providing  a  model  for  the  Total 
Ship  System  Engineering  (TSSE)  program  at  the  Naval  Postgraduate  School.  The  TSSE 
program  produces,  as  a  capstone  design  project,  one  proposed  ship  each  year  by  the  students 
in  the  program.  The  ships  designed  by  the  TSSE  program  currently  have  no  maneuverability 
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prediction  program.  It  is  hoped  that  the  final  product  of  this  study  will  provide  the  TSSE 
design  teams  with  a  tool  enabling  further  study  of  their  proposed  ships  form  the  standpoint 
of  maneuverability  and  station-keeping. 

B.  OBJECTIVES 

The  goal  of  this  study  is  to  construct  a  ship-steering  simulation  model,  complete  with 
power  plant  dynamics  and  wind  and  wave  disturbances,  which  will  be  used  to  assess 
alternative  control  laws  for  autonomous  station-keeping  during  at-sea  replenishment.  This 
simulation  will  ultimately  be  provided  to  the  TSSE  program  for  use  by  future  ship  design 
teams.  Specific  objectives  include: 

1) .  The  development  of  a  mathematical  model  for  ship  dynamics,  power  plant  and 
rudder  dynamics,  and  disturbance  factors  for  wind  and  wave  conditions  in  various  sea  states. 
The  relevant  mathematical  equations  appear  in  Chapter  II  below. 

2) .  The  development  of  a  computer  simulation  based  on  the  equations  derived  in  the 
mathematical  model.  This  simulation  development  is  described  in  Chapter  HI  below. 

3) .  The  development  of  a  controller  design,  as  described  in  Chapter  IV  below. 

4) .  Test  and  evaluation  of  one  controller  design  using  the  simulation.  These  results 
are  discussed  in  Chapter  V  below. 

Input  for  the  simulation  study  will  be  taken  from  the  final  design  report  for  the  first 
ship  designed  by  the  TSSE  program:  RDS-2010. 
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H.  MATH  MODEL  DEVELOPMENT 


A.  COORDINATE  SYSTEMS 

Two  coordinate  systems  are  used  in  the  description  of  the  ship's  movement  relative 
to  the  earth.  The  first  is  a  system  fixed  with  respect  to  the  earth  oriented  at  some  fixed  point 
on  the  surface;  the  other  is  a  system  fixed  with  respect  to  the  three  major  axes  of  the  moving 
ship.  A  ship-fixed  coordinate  system  has  the  advantage  of  constant  moments  of  inertia  and 
constant  moment  arm  lengths  in  all  three  directions.  An  earth-fixed  coordinate  system,  if 
used,  would  constantly  change  the  lengths  of  moment  arms  in  each  direction  as  the  ship  turns 
and  rolls.  Convention  applies  control  and  disturbance  forces  to  the  ship  in  ship-fixed 
coordinates,  resulting  in  motion  of  the  ship's  coordinate  origin  relative  to  the  earth.  The 
motion  of  the  coordinate  axis  of  the  ship  represents  the  ship's  movement  with  respect  to  the 
earth.  Figure  1  illustrates  the  coordinate  system  being  used.  It  is  a  right  hand  system 
positive  forward,  starboard,  and  down,  using  x,  y,  and  z,  respectively,  to  denote  the  ship 
fixed  axes.  Also  shown  are  the  Xo  and  Yo  axes,  which  are  considered  fixed  in  the  earth. 
Xo(t)  and  Yo(t)  display  a  distance  the  ship  has  traveled  with  respect  to  the  earth  origin  point 
at  a  time  (t).  The  rudder  angle  8  is  shown  as  convention  employs;  positive  rudder  is  for  a  port 
turn.  The  angle  *¥  represents  the  angle  the  ship's  velocity  vector  is  displaced  from  a  parallel 
of  the  Xo  axis  placed  at  the  ship's  center  of  gravity.  For  simplicity,  the  ship's  velocity  vector 
will  be  considered  to  be  along  its  longitudinal  axis,  making  the  angle  representative  of  the 
ship's  heading  with  respect  to  the  earth  fixed  origin  (X^Y^). 
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Figure  1.  Earth  Centered  and  Ship  Centered  Coordinate  Systems 

B.  EQUATIONS  OF  MOTION 


An  abbreviated  presentation  is  made  of  a  ship's  equations  of  motion.  A  more  in  depth 
discussion  of  these  equations  can  be  found  in  Principles  of  Naval  Architecture. 

[Ref.  2:Ch.  8] 

1.  Earth  Fixed  Equations 

From  Principles  of  Naval  Architecture  [Ref  2],  a  ship's  movement  with  respect 
to  an  earth  fixed  reference  point  (X^Y^)  can  be  described  by  equations  (la  and  b)  as: 

(a)  X0G  =  ucosy  -  vsiny 

( b )  Yog  =  vcos\j/  +  wsinvj/ 
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where:  x,  y,  and  z  are  in  ship  fixed  coordinates,  and  u  =  x,  v  =  y,  \j/  =  ship 's  heading, 
are  termed  surge,  sway  and  yaw  respectively. 

2.  Ship  Fixed  Equations 

If  a  moving  set  of  axes  is  desired,  Ref  2  shows  that  a  ship's  equations  of 
motion  on  a  horizontal  plane  are: 

(a)  Xu  (w  -  ux)  +  (m  -  X.)  u  =  0 

(b)  -Yvv  +  (m-Y.)v-  (Yr  -  m)  r  -  (Yr  -  mxG)  f  =  Ys 5  (2) 

(c)  -Nv  -  (N.  -  mxG)  v  -  (Nr  -  mxQ)  r  +  (Jz  -  N.)  f  =  Ns5 

where:  X  =  total  force  in  the  x  direction 
Y  =  total  force  in  the  y  direction 
N  =  total  moment  about  the  z  axis 
I2  =  mass  moment  of  inertia  about  the  z  axis 
m  =  ship's  mass 

Xg  =  distance  ship's  center  of  gravity  is  displaced 
from  the  centerline  (positive  forward) 

8  =  rudder  deflection  angle 

r  =  'F 

Equations  (2a-c)  assume  symmetry  about  the  ship's  longitudinal  axis,  and  employs 
a  shorthand  notation  where: 
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(a)  Yv  =  8Y/dv 

(b)  N.  =  8N/dr,  etc.... 


(3) 


Assuming  linearity  of  the  above,  we  see  that  X  in  equation  (2a)  is  not  coupled  with 
second  two  equations  (2a  and  b),  hence  it  will  be  dealt  with  separately. 

a.  Nondimensionalizing  Variables 

Equations  of  (2a  and  b)  has  units  offeree,  while  equation  (2c)  has  units  of  a 
moment.  Therefore,  it  is  convenient  to  nondimensionalize  the  variables  as  follows: 


( a ) 
(*) 
(P) 

(d) 

(e) 
if) 
is) 
ih) 
(0 
(/•) 


_  /  _ 
m  = 


m 


p/2  V 
jf  m 

2  p/2  L5 

\7 

yL  = 


p/2  L2U2 


Y'  = 


Yi  = 


Y !  = 


n!  = 


n!  = 


K  = 


n!  = 


p/2  L3U 
p/2  L3U 

h 

p/2  L4U 

K 

p/2  L3U2 

K 

p/2  Z4/7 

p/2  X4t/ 
N. 


p/2  L 


(4) 
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where:  p  =  density  of  seawater 


L  =  ship's  length 
U  =  ship  operating  speed 

b.  Formulae  for  Predicting  the  Hydrodynamic  Derivatives 
From  Regional  Deterrence  Ship  (RDS-2010)  [Ref.  l:Ch.  8:  Sec.  10],  the 
following  formulae  are  found  for  deriving  hydrodynamic  derivatives  from  hull  data: 

(«)  <  --  wl\  *  (Xl), 

(b)  1/2  (7vV 

w  ri  -  <rX  +  pIvot), 

(d)  7,  -  0  -  1/2(74 

(«)  i  mxi>, 

(0  If!  =  <K\  *  1/4 (7„V  (S) 

(?)  a//  =  +  mriood, 

(h)  n!  =  0  -  1/2(74 

©  4,  =  (i4 

CO  JVjtf  -  -1/2(74 

where:  ( )h  denotes  the  contribution  to  the  variable  from  the  hull,  and  ( )  f  denotes  the  sum 
of  contributions  to  the  variable  from  fins  and  other  appendages  as  shown  in  equations  (6a-g). 
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where: 


(a)  (Y\  =  -UTIL  +  CD 

(b)  (Y')h  =  kxm  ’  +  xp/L(Y\ 

ip)  +  xJUJX 

(d)  QfX  =  x'L  +  (#4 

(*)  vX  =  r  «o 

bow 

stem 

W  ^  =  *7171  /  c/*v* 

^  *  bow 


(«) 

(*) 

(0 


stern 


m'z  =  (h'/k^)  m2 


stem 

f  Csd2  xdx 

bow 

stem 

/  c^Ia 

bow 


(7) 
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and  as  empirically  derived  by  Vann  [Ref.  3]: 


(а)  kx  =  0.3(2  TIL) 

(б)  k2  =  1.0  -  0.5(2 TIL) 

(c)  k'  =  1  -  1.333(277L) 

id)  Cs  =  CQ  +  Cx(sa/bd)  +  C2(salbd)2  (8) 

(e)  C0  =  0.8572  +  0.5330(4^/6) 

(/)  C,  =  3.374  -  1.3661(4^6) 

(g)  C2  =  -1.7323  +  Q.Z61Q(4dld) 

where: 

T=  ship's  draft 

CD=  ship's  coefficient  of  drag 
kj  =  longitudinal  coefficient  of  accession  to  inertia 
k2  =  lateral  coefficient  of  accession  to  inertia 
k'  =  rotational  coefficient  of  accession  to  inertia 
Cs  =  two  dimensional  sectional  inertia  coefficient,  calculated  using  strip 
integration  along  the  ship’s  hull, 
sa  =  section  area  of  section  being  considered 
d  =  draft  at  section  being  considered 
b  =  beam  at  section  being  considered 
Co ,  C„  C2  =  interim  variables  used  in  the  calculation  of  Cs 
Xp  =  distance  from  centerline  to  point  of  application  of  fluid  force 
x^L  =  14  prismatic  coefficient 
Af  =  profile  area  of  the  appendage 
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Af'  =  nondimensionalized  profile  area  of  fin  or  appendage  (A/LT) 
a  =  aspect  ratio  of  the  fin  or  appendage,  defined  by  Gillmer  [Ref.  4]  as: 


h  =  span  of  the  appendage 

Calculation  of  the  above  derivatives  is  accomplished  in  MATLAB  code  with  the  files 
HYDROGEN.M  and  CSGEN.M  that  are  listed  in  appendices  A  and  B. 

3.  Propulsion  Equations 

To  solve  for  motion  along  the  ship's  longitudinal  axis,  equations  (10a  and  b)  for  gas 
turbine  engines  are  used,  as  presented  by  Tozzi.  [Ref.  5] 

(a)  V  =  g/m(Tp  -  Rs) 

(i)  *  =  -Ug,  -  Qf  -  gp  (10) 


where:  V  =  ship  velocity  in  feet  per  second 

N  =  Propeller  shaft  speed  in  rotations  per  second 

g  =  acceleration  due  to  gravity  in  feet  per  second  squared 

m  =  ship  mass  in  pounds  mass 

Tp  =  propeller  thrust  in  pounds  force 

Rs=  resistance  of  the  ship  in  pounds  force 


I  =  moment  of  inertia  of  the  drive  train  referred  to  the  propeller  shaft 
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Qt  =  engine  torque  available  at  the  propeller  shaft  in  foot  pounds  force 
Qf  =  friction  torque  losses  in  foot  pounds  force 
Qp=  propeller  torque  losses  in  foot  pounds  force 

a.  Propeller  Thrust 

The  thrust  is  a  measure  of  force  developed  by  the  propeller  in  the  direction  of 
the  shaft.  Mathematically,  it  is  expressed  as: 

TP  =  kT  p  D4  N 2  (11) 


where:  kx  =  thrust  coefficient 

p  =  density  of  seawater  in  slugs  per  cubic  foot 
D  =  propeller  diameter  in  feet 


b.  Ship  Resistance 

The  resistance  a  hull  form  has  is  normally  found  during  scale  model  tests  in  a 
wave  tank.  For  the  RDS-2010  the  resistance  is  estimated  as  a  function  of  speed.  From  data 
provided  by  Alexander  [Ref.  1]  it  is  seen  that  ship  resistance  varies  as: 


Rs  =  172  V2 


(12) 


c.  Drive  Train  Moment  of  Inertia 

As  presented  by  Tozzi  [Ref.  5],  the  moment  of  inertia  includes  the 
contributions  of  all  drive  train  components  referred  to  the  propeller  shaft.  For  this  simulation. 


11 


Equation  (13)  is  used  to  estimate  the  moment  of  inertia,  based  on  Ref.  5. 


/  =  2.0  x  105 


(13) 


d.  Shaft  Torque 

The  shaft  torque  (Qx)  as  supplied  in  the  report  on  the  RDS-2010  includes  the 
effect  of  the  reduction  gears  in  the  drive  train,  and  is  a  measure  of  the  torque,  in  foot  pounds 
force,  provided  to  the  propeller  shaft. 

e.  Friction  Torque 

Friction  in  the  entire  drive  train  can  be  modeled  as  a  loss  in  torque  developed 
to  produce  thrust.  From  Ref  2  it  is  seen  that  friction  torque  can  be  approximated  as: 

Qf  =  6000  N  (14) 


f.  Propeller  Torque 

The  torque  required  to  rotate  a  propeller  under  various  conditions  is  normally 
measured  during  open  water  tests.  For  the  subject  ship  it  can  be  calculated  using  data 
provided  in  Ref.  1  as: 

Qp=Kq  p  N2  Ds  (15) 

where:  Kq  =  Torque  coefficient  for  the  hull  under  consideration. 
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4.  Steering  Equations 

A  low  order  rudder  model  can  be  developed  similar  to  Ref.  6.  It  is  modeled  as  an 
integrator  with  gain  Kg'  as  follows: 


(«) 

(*) 


u 


K  = 


(16) 


where:  5m  =  maximum  rudder  angle 
8einax=  maximum  error  input. 


5.  Ship  Dynamics 

To  convert  the  Y  and  N  equations  of  Equations  (2a  and  b)  to  a  form  usable  in  the 
simulation  it  is  necessary  to  develop  matrices  as  follows: 


(a)  [M]  x  =  [H]  x  +  [R]  u 

(b)  x  =  [M]'1  [H]  x  +  [My1  [i?]  u 


(17) 
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where: 


(«) 

(b) 

ip) 

id) 


[Ml 

[*] 


(»'  -  -  *”/) 
(«'*,'  -  <x/;  -  jv/) 

[  yXy;  -  wo 

[jNT'CAT/  m^) 

>/] 


x  =  [v  r]r 
u  =  8 


(18) 


The  new  matrices  formed  are  implemented  as  follows: 


id) 

x  =  Ax  +  Bu 

ib) 

g 

*— t 

5 

u_Li 

11 

(19) 

ic) 

B  =  lM]-‘  [tf] 

C.  DISTURBANCES 

To  test  the  system's  performance  under  varying  conditions,  wind  and  sea  disturbances 
are  added  to  model  the  effect  of  increasing  sea  state  on  the  system.  Also  included  is  the 
venturi  effect,  a  phenomenon  affecting  the  sway  forces  and  yaw  moments  of  ships  as  they  pass 
each  other  in  close  quarters. 

Lastly,  measurement  noise  is  added  to  the  observation  of  the  supply  ship's  position  to 
model  the  effect  of  winds  and  seas  acting  on  the  supply  ship.  This  noise  models  the  effect  of 
sea  state  by  increasing  position  estimate  error  as  sea  state  increases. 
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1.  Sea  Disturbances 


In  order  to  simply  model  what  can  be  made  infinitely  complex,  Uhrin  [Ref.  7]  employs 
formulas  (20a-f)  for  a  wave.  These  formulas  generate  a  disturbance  at  a  primary  frequency 
and  its  second  harmonic: 

(a)  W  =  WF  (1  +  WRV)  sin  (WEF)  +  sin  (2  WEF) 

(b)  mF  =  ink (u  +  ws  cos  (XPrs) 

(c)  -  WSff/p  (20) 

(d)  Xs  -  W  costae 

(e)  Ys  =  W  sinT^ 

if)  Ns  =  W  sin  (2 

where:  W  =  total  wave  force 

WEF  =  wave  encounter  frequency 

WF  =  maximum  wave  force 

WL  =  wave  length,  normalized  by  WLT/LOA 

WLT  =  true  wave  length 

LOA  =  ship  length  overall 

WS  =  wave  speed  in  feet  per  second 

Trs  =  relative  seas  heading 

Tts  =  true  seas  heading 

'J'shp  =  ship  heading 

X,,  =  sea  force  in  x  direction 
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Ys=  sea  force  in  y  direction 
Ns  =  sea  moment  in  N  direction 
WRV  =  wave  random  variable 

The  wave  random  variable  is  added  to  model  the  randomness  of  the  sea.  It  is 
implemented  as  zero  mean  gaussian  white  noise  with  a  variance  of  0.01. 


2.  Wind  Disturbances 

To  account  for  the  wind's  effect  on  the  ship,  a  model  is  designed  similar  to  that 
employed  by  Clark.  [Ref.  8]  It  is  necessary  to  take  wind  true  velocity  and  speed  in  the  earth 
fixed  reference  frame,  and  translate  each  component  of  the  forces  and  moment  to  the  ship 
fixed  reference  frame. 

First,  the  relative  velocity  and  direction  must  be  found: 

(a)  Vrw 

0)  yRiy 

Vw  cos  (yw  VsHip)  u 


Ml  ~ 


(21) 


Coefficients  for  each  component  of  wind  effect  are  developed  as  follows: 


(«) 

0) 

(c) 


PaAf 

7000pw  LBP 2 
-  Pajs 
8000pw  LBP 2 
Pa  4  LOA 
Pw  LBP 3 


(22) 
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Lastly,  the  nondimensional  forces  and  moment  are  developed  in  equations  (23a-c): 

(a)  Xlr  =  (-^)2  sin[9/7(|Vw  |  -  n)] 

(*)  sin(¥^) 

(c)  sin(2^) 

where:  Xw'  =  wind  force  in  x  direction 
Yw’  =  wind  force  in  y  direction 
Nw'  =  wind  moment  about  yaw  axis 
Yrw  =  relative  wind  velocity 
^rw  =  relative  wind  direction 
Vw  =  true  wind  velocity 
VPTW  =  true  wind  direction 
'f'sHjp  =  ship  heading 
u  =  ship  surge 
v  =  ship  sway 

Cwx  =  x  direction  wind  coefficient 
Cwy  =  y  direction  wind  coefficient 
Cwn  =  N  direction  wind  moment 
pa  =  mass  density  of  air 
pw  =  mass  density  of  water 


(23) 
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Af  =  ship  frontal  area 

A*  =  ship  side  area  (estimated  as  10Af ) 

LO  A  =  ship  length  overall 

LBP  =  ship  length  between  perpendiculars 

Xw'  =  nondimensional  x  direction  wind  force  component 

Yw'  =  nondimensional  y  direction  wind  force  component 

Nw‘  =  nondimensional  N  direction  wind  force  component. 

3.  Venturi  Disturbances 

The  venturi  effect  can  be  closely  estimated  as  a  nonlinear  function  of 
longitudinal  separation  that  is  inversely  proportional  to  lateral  separation  as  seen  in  Ref.  6. 
Equations  (24a  and  b)  are  used  to  approximate  this: 

(a)  y!^  =  0.45)  «  10-5  [smc  (0.0029  IS)]  Q 

„/  _  _  „/  <24) 

\°)  iyVENT  ~  1  VENT 


where: 

Yvent'  =  nondimensional  venturi  force 
Nvent'  =  nondimensional  venturi  moment 
LS  =  longitudinal  separation  in  feet 

Q  =  lateral  separation  multiplier  ( varies  linearly  from  1  to  0.5  as  separation  increases 
from  100  to  200  feet). 
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4.  Measurement  Noise 


The  measurement  noise,  as  described  above,  is  simply  modeled  as  gaussian 
white  noise  with  a  variance  of  5,  10  or  20,  which  will  vary  with  sea  states  0,  3,  and  6, 
respectively. 
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m.  SIMULATION  DEVELOPMENT 


A.  PROBLEM  GEOMETRY 

The  underway  replenishment  scenario  considered  in  this  research  will  place  the  supply 
ship  at  an  earth  fixed  origin  (Xo,Yo)  at  problem  start.  Throughout  the  problem  the  supply 
ship  will  remain  at  the  relative  origin.  The  receive  (subject)  ship  will  start  the  problem  in 
waiting  station,  one  thousand  yards  astern  of  the  supply  ship.  Both  ships  will  begin  the 
problem  on  base  course,  at  base  speed  (course  000  at  14  knots).  This  system  will  drive  to  a 
starboard/  port  supply/  receive  geometry.  The  system  will  seek  to  drive  the  receive  ship  to 
a  position  one  hundred  and  thirty  feet  abreast  of  the  supply  ship,  at  (0,130)  relative  position. 


Ysd  =  y  separation  desired  = 130 

_  _ i 

I 

1 - 1  Yr(t)  =  130-yr(t):  y  remaining 

S(t) 

. 4(0 

\^tj )  D(t)  =  sqrt[xr(t)*2  +  (130-yr(t))*2] 

So 

/  Ao 

s'  Do  =  sqrtf 3000*2  +  130*2] 

Vy>^R(t)  =  (xr(t),yr(t)) 

Ro  =  (-3000,0) 

Figure  2.  Problem  geometry. 


The  initial  positions  of  the  supply  and  receive  ships  is  displayed  in  Figure  2  as  S0  and 
Ro  respectively.  The  initial  aim  point  is  denoted  as  A„.  From  the  figure  it  is  evident  that  the 
aim  point  moves  with  the  supply  ship  at  base  course  and  speed.  As  the  problem  progresses. 
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the  ship  is  at  a  relative  position  (xr(t),yr(t)),  steers  an  angle  v|/(t),  and  has  a  distance  to  travel 
D(t).  Thus,  the  initial  angle  vj/0  and  distance  D0  are  not  constantly  fixed  with  respect  to  the 
absolute  origin,  and  must  be  computed  and  fed  into  the  system  as  command  signals  for  all 
time  (t). 

B.  SHIP'S  PLANTS 

1.  Rudder 

As  previously  mentioned,  a  low  order  model  for  a  rudder  plant  was  designed  in  Ref. 
6,  and  is  shown  as  Figure  3.  Heading  error  is  input  to  the  system  at  the  in  block.  The  rudder 
stops  limit  the  amount  of  travel  and  is  modeled  as  a  saturation  block.  The  rate  limit  block 
simulates  a  hydraulic  pump  that  limits  the  amount  of  error  into  the  rudder,  which  is  modeled 
as  an  integrator  with  gain  Kg’.  The  rudder  angle  (5)  is  the  output  of  the  subsystem  at  the  out 
block.  The  SIMULINK  implementation  requires  renaming  the  variable  Kg'  in  the  rudder 
subsystem  as  Kgpri. 


Figure  3.  Rudder  subsystem. 
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2.  Ship's  Engine 

The  Engine  model  is  designed  to  output  the  thrust  necessary  to  drive  the  ship.  Input 
is  the  desired  torque.  Next,  the  subsystem  develops  the  shaft  acceleration  as  defined  earlier. 
An  integrator  is  used  to  generate  the  propeller  rotations  per  second,  which  is  then  converted 
to  a  thrust  value.  The  initial  condition  set  on  the  integrator  starts  the  problem  with  the  ship 
already  at  base  speed.  The  thrust  output  is  then  calculated  as  a  function  of  shaft  rotations 
as  noted  previously.  From  Equations  (10a  and  b)  the  conversion  into  a  SIMULINK 
environment  is  relatively  straightforward  as  seen  in  Figure  4. 


Figure  4.  Engine  Subsystem. 
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C.  SHIP’S  DYNAMICS 


1.  Steering  Dynamics 

To  model  the  ship  dynamics  in  SIMULINK,  the  matrices  in  Equations  (19a-c)  are 
used  in  the  subsystem.  These  matrices  lend  themselves  very  well  to  a  SIMULINK 
implementation.  The  subsystem  also  adds  the  Y  and  N  disturbances  to  the  accelerations. 
Figure  5  displays  the  ship  dynamics  subsystem. 


Figure  5.  Steering  dynamics  subsystem. 


2.  Propulsion  Dynamics 


Equation  (10A)  is  used  to  produce  the  ship's  longitudinal  velocity  from  the  developed 


propeller  thrust.  Again,  the  conversion  is  uncomplicated  as  seen  in  Figure  6. 
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Figure  6.  Propulsion  dynamics  subsystem. 


D.  EARTH  ORIENTED  POSITION  GENERATION 


Equations  (la  and  b)  are  very  simply  translated  into  a  SIMULINK  subsystem  to  solve 
for  velocities  in  the  Xo  and  Yo  directions.  These  velocities  are  integrated  to  obtain  the 
positions  with  respect  to  the  earth  fixed  origin.  Figure  7  shows  the  earth  position  generation 
subsystem  as  implemented. 
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Figure  7.  Earth  oriented  position  generation  subsystem. 


E.  RELATIVE  POSITION  GENERATION 

In  order  to  generate  the  ship's  position  relative  to  the  supply  ship,  the  position  of  the 
supply  ship  with  respect  to  the  earth  fixed  origin  must  be  compared  with  the  receive  ship's 
position  with  respect  to  the  same  origin  at  each  time  interval.  This  is  accomplished  by  using 
base  speed  as  the  ideal  supply  ship  motion  along  the  Xo  axis.  The  supply  ship's  motion  along 
the  Yo  axis  is  ideally  zero.  To  both  position  estimates  of  the  supply  ship,  the  measurement 
noise  is  added.  This  effectively  models  a  ship's  sensor  (radar),  and  accounts  for  the  supply 
ship's  motion  in  an  increasing  sea  state.  The  relative  position  generator  subsystem  is  shown 
in  Figure  8  as  implemented. 
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Figure  8.  Relative  position  generation  subsystem. 


F.  DISTURBANCE  GENERATION 


As  discussed  in  the  previous  chapter  the  disturbances  are  to  test  the  system's 
performance  under  less  than  ideal  conditions.  The  forces  and  moments  generated  are  summed 
from  each  component. 
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1.  Wind  Disturbance 


From  the  previous  chapter  the  Figure  9  shows  the  subsystem  developed  to  generate 
wind  disturbances. 


Figure  9.  Wind  disturbance  subsystem. 
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2.  Sea  Disturbance 


The  mathematical  model  for  seas  as  developed  previously  translates  very  well  into  a 
SIMULINK  subsystem,  and  is  shown  in  Figure  10. 
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3.  Venturi  Disturbance 


The  venturi  phenomenon  is  modeled  in  SIMULINK  with  the  subsystem  shown  as 


Figure  1 1 .  The  subsystem  follows  directly  from  Equations  (27a  and  b). 


Figure  11.  Venturi  disturbance  subsystem. 
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IV.  CONTROLLER  DESIGNS 


Two  types  of  controllers  are  employed  within  the  system  designed.  The  first  is  a 
command  guidance  signal  styled  to  drive  the  ship  similar  to  a  conning  officer  during  an 
underway  replenishment.  This  is  accomplished  by  the  conversion  of  the  receive  ship's  position 
relative  to  the  supply  ship  into  command  signals  to  the  rudder  and  engine.  The  second  is  a 
root  locus  design  for  the  engine  plant  to  drive  the  ship  more  efficiently. 

A.  COMMAND  SIGNAL  GENERATION 

This  first  controller  design  is  implemented  as  a  separate  subsystem,  converting  the 
receive  ship's  relative  X  and  Y  positions  into  rudder  and  engine  command  signals  as 
mentioned  above.  In  order  to  achieve  the  desired  goal  of  imitating  a  conning  officer,  the  first 
goal  is  to  have  the  ship  open  the  lateral  separation  in  order  to  safely  close  the  longitudinal 
distance  without  danger  of  collision.  While  the  rudder  steers  the  ship  towards  the  aim  point, 
the  engine  must  be  sped  up  to  close  the  longitudinal  separation.  A  typical  conning  officer  will 
bring  the  ship  up  to  a  speed  of  20  knots  to  expeditiously  gain  station.  Once  the  ship  is  within 
approximately  200  feet  of  station,  the  ship  is  slowed  to  1  knot  above  base  speed  to  avoid 
going  past  the  station. 

Once  the  ship  is  in  station,  the  goal  is  to  stay  there.  This  requires  small,  continual 
adjustments  in  both  the  rudder  and  the  engine  in  response  to  the  effects  of  the  wind,  sea,  and 
the  presence  of  the  supply  ship. 

To  gain  and  then  maintain  station  defines  two  separate  problems  as  defined  above. 
The  first  is  an  approach  problem  requiring  large  corrections.  The  second  is  a  station-keeping 
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problem  requiring  finer  adjustments  to  maintain  the  desired  position  and  avoid  dangerous 
maneuvers  while  in  close  quarters  to  the  supply  ship. 

During  the  approach  phase,  the  desired  heading  is  conveyed  to  the  rudder  subsystem 
as  'F(t),  as  defined  previously.  The  distance  D(t)  is  used  to  generate  a  velocity  command 
signal  ordering  the  engine  to  20  knots  until  within  200  feet  of  station.  The  velocity  command 
signal  is  also  converted  to  a  torque  command  signal  based  on  the  predicted  torques  at  various 
speeds  for  the  subject  ship. 

The  lateral  separation,  being  opened  to  the  proper  distance,  is  employed  as  the 
division  between  approach  and  station-keeping  phases.  Once  the  lateral  distance  is  within  5 
feet  of  the  desired  setting,  the  command  generation  subsystem  toggles  to  the  station-keeping 
phase.  During  this  phase  the  rudder  is  restricted  to  a  maximum  amplitude  of  one-fifth  of  'F0. 
The  rudder  will  steer  the  ship  with  minute  adjustments  when  within  20  feet  of  the  desired 
lateral  separation.  The  station  keeping  phase  also  restricts  the  engine  response  by  basing  the 
D(t)  on  only  the  longitudinal  separation.  Further,  the  velocity  command  is  regulated  to  within 
1  knot  of  base  speed  when  within  40  feet  of  station.  The  engine  will  respond  to  gradually 
speed  up  or  slow  down  in  response  to  the  relative  longitudinal  position.  The  rather 
complicated  verbal  description  above  is  actually  very  simple  to  implement  in  SIMULINK. 
The  relative  positions  are  used  as  input  and  converted  to  'f'(t)  and  D(t)  signals  as  described 
in  the  problem  geometry.  These  are  used  with  look-up  tables  to  accomplish  the  intent  of  the 
subsystem.  The  subsystem  is  pictured  in  Figure  12  as  implemented. 
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Figure  12.  Command  signal  generation  subsystem. 


B.  ROOT  LOCUS  CONTROLLER 

The  engine  requires  a  separate  controller  developed  to  improve  the  response  of  the  plant. 
Figure  13  displays  the  actual  plant  implemented. 
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Y  =  u 


Figure  13.  Actual  engine  plant. 
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where: 


<«) 


S  X  = 


u  -  uQ 
N-N0 


(P)  5  Qt^Qt-  QTo 


(27) 


The  first  term  in  the  above  non-linear  Equation  (27b)  is  assumed  to  go  to  zero,  leaving 
a  linear  equation.  The  second  term  defines  the  plant's  A  matrix.  The  third  term  defines  the 
plant's  B  matrix.  For  the  controller  design,  a  new  5QT  is  defined  as  below. 

5 Qt  =  K  (u  -  w0)  =  K  Sv  (28) 


Then  for  the  approximate  system: 

5x  =  A  8x  +  B  K  5u  =  A  fix  +  B  K  C  fix  (29) 


hence: 


(a)  5x  =  [A  -  KB  C\  hx 
(i b )  5y  =  C  5x 

because  K  is  a  scalar. 

The  poles  of  this  system  are  at: 

DET(sI  -  A  +  KBC ) 


(30) 


(31) 


and,  the  system  transfer  function  is: 

C  [si  -  A  =  +  KBC]' 1  B 


(32) 
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A  suitable  value  for  K  was  found,  experimentally  to  be  100000.  This  was  employed 


in  the  system  as  displayed  in  Figure  14. 


Desired 

Velocity 

Input 


Figure  14.  Engine  plant  as  implemented. 
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V.  SIMULATION  RESULTS 

The  system  was  run  several  times  to  demonstrate  its  performance  under  varying 
conditions.  First  the  system  was  run  without  disturbances  as  shown  in  Figures  15-17..  Next, 
the  system  was  run  with  only  the  measurement  noise  applied  as  a  disturbance  seen  in  Figures 
18-20.  For  this  run  the  measurement  noise  was  increased  to  20  to  better  display  its  effects 
on  the  system.  Lastly,  the  system  was  run  with  wind,  sea,  and  venturi  disturbances  applied 
in  addition  to  the  measurement  noise  shown  in  Figures  20-23.  For  the  last  run  the  wind  and 
seas  were  set  at  5  knots,  the  normalized  wavelength  set  at  0.1,  and  the  true  direction  set  at 
000  degrees  true  (on  the  bow).  For  each  run  three  plots  are  shown.  The  first  plot  presents 
the  relative  trajectory  of  the  receive  ship  (a  bird's  eye  view  of  the  track  taken  by  the  receive 
ship).  The  next  two  plots  display  the  lateral  and  longitudinal  trajectories  plotted  against  time 
to  compare  settling  times  of  the  different  cases. 


Figure  15.  Undisturbed  system  relative  trajectory 
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Figure  17.  Undisturbed  system  longitudinal  trajectory. 
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Figure  19.  Measurement  disturbance  lateral  trajectory. 
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Figure  20.  Measurement  disturbance  longitudinal  trajectory. 
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Figure  22.  Full  disturbances  lateral  trajectory. 
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VL  CONCLUSIONS  AND  RECOMMENDATIONS 


A.  CONCLUSIONS 

From  the  simulation  results  in  the  previous  chapter  it  is  evident  that  the  system 
performs  as  desired  when  undisturbed.  The  lateral  separation  is  opened  safely  while  closing 
the  longitudinal  distance.  And,  once  the  station  is  attained,  the  system  keeps  the  ship  at  the 
desired  point.  Additionally,  it  is  shown  that  the  measurement  noise  has  little  effect  on  the 
system.  In  each  of  these  first  two  cases  the  system  has  a  time  to  station  of  420  seconds. 

However,  it  is  also  evident  that  the  system  does  not  comfortably  open  the  lateral 
separation  while  being  exposed  to  the  wind,  sea,  and  venturi  disturbances.  But,  the  system 
still  has  a  satisfactory  time  to  station  of  430  seconds,  which  is  very  comparable  to  the 
undisturbed  case.  And,  once  on  station  the  system  still  keeps  the  ship  in  place. 

B.  RECOMMENDATIONS  FOR  FURTHER  STUDY 

Based  on  the  above  conclusions,  the  first  recommendation  is  to  design  a  controller 
to  better  handle  the  rudder  during  exposure  to  disturbances.  Further,  many  different  control 
schemes  for  both  the  engine  and  the  rudder  may  be  designed  and  compared  using  classical 
as  well  as  optimal  control  approaches. 

As  the  TSSE  program  has  designed  two  additional  ships  since  the  RDS-2010,  it  may 
be  desired  to  apply  this  study  to  the  other  ships.  This  may  require  the  development  of 
different  propulsion  plants  which  could  be  stored  and  modified  for  future  use  by  the  TSSE 
program.  Additionally,  the  TSSE  program  could  employ  this  system  to  test  more 
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fundamental  performance  measures  of  their  proposed  ships,  such  as  tactical  diameter.  It 
should  also  be  noted  that  this  system  would  be  applicable  to  a  commercial  vessel  design. 

Finally,  the  disturbance  model  would  prove  useful  to  any  design  team  wishing  to  test 
the  sturdiness  of  a  proposed  hull  design.  The  forces  and  moments  the  hull  would  be  exposed 
to  can  be  measured,  and  compared  to  the  maximum  levels  the  structure  should  be  exposed 
to. 
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APPENDIX  A  HYDROGEN.M 


%  HYDROGEN .M:  FILE  TO  GENERATE  HYDRODYNAMIC  DERIVATIVES,  WIND 
%  COEFFICIENTS ,  WIND  AND  SEAS  TRUE  DIRECTION,  AND  ENGINE 
%  SIMULATION  CONSTANTS  FOR  SIMULINK  FILE  "SHIPSIM". 

%  VALUES  REQUIRING  USER  INPUT  FOR  SYSTEM  CUSTOMIZATION  ARE 
%  NOTED  AS  (input)  IN  THE  COMMENT  LINE. 

%  THIS  FILE  MUST  BE  RUN  PRIOR  TO  "SHIPSIM"  TO  GENERATE  THE 
%  NECESSARY  INPUTS 


rho= 1 .99 14*0.03 11;%  Density  of  H20  in  Slugs/Cub.  Ft. 

%  at  55  deg  F  *  (Slugs/Lb  Mass) 

U  =  14*  1 .689  ;%  Nominal  Ship's  FWD  Speed 

%  (14  knots)  in  Feet  per  Sec. 

%  HULL  COEFFICIENT  OF  DRAG  AT  ZERO  ANGLE  OF  ATTACK  CALCULATION 

rhos=  1 .99 14  ;%  Density  of  H20 

vs= 1 .3535  ;%  Kinematic  Viscosity  of  H20 

q=.5*rhos*vsA2  ;%  Dynamic  Pressure  of  H20 
Do=123559  ;%  Drag  Force  at  16  Knots 
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WSA=24204. 1  ;%  Wetted  Surface  Area  at  Displacement 

Cdo=Do/(q*WSA)  ;%  Hull  Coefficient  of  Drag 

%  NECESSARY  INPUTS  FOR  HYDRODYNAMIC  DERIVATIVE  GENERATOR 


D  =5721.8 

;%  (input)  Ship  Disp.  in  Long  Tons  (5721.8) 

T  =  15.01 

;%  (input)  Ship's  Draft  in  Feet  (15.01) 

L  =409.315 

;%  (input)  Ship’s  Length  (LOA)  (409.315) 

B  =  55.09 

;%  (input)  Ship's  Beam  in  Feet  (55.09) 

Cb  =0.599 

;%  (input)  Ship's  Block  Coefficient  (0.599) 

mship  =  (2240*D)/32.2  ;%  Ship's  Mass  Calculator 

delmaxdot  =  2 

;%  (input)  maximum  rudder  rate 

delemax  =  7 

;%  (input)  maximum  error  input 

Kg  =  delmaxdot/delemax  ;%  Rudder  Gain  Calculator 

Xg  =-5.85 

;%  (input)  Estimate  of  Xg  (LCB  for  Disp.) 

xp  =-22.5 

;%  (input)  Long.  Center  of  Float.  (LCF) 

%  at  Disp. 

Priz=  0.650 

;%  (input)Prismatic  Coefficient 

%  APPENDAGE  DATA:  CALCULATED  FOR  RUDDERS  AND  THE  SONAR  DOME, 
%  MAY  BE  CUSTOMIZED  FOR  OTHER  APPENDAGES  BY  ADDING  THE 
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NECESSARY 


%  CALCULATIONS  SIMILAR  TO  THOSE  BELOW. 

ProfArudl  =  1 17.4  ;%  (input)  Rudder  Profile  Area  (1 17.4) 

PAnondl  =  ProfArudl/(L?T)  ;%  Non  Dimensionalized  Profile  Area 

hrudl  =11.95  ;%  (input)  Rudder  Height 

arudl  =  hrudl A2/ProfArudl  ;%  (input)  Rud  Eff  Aspect  Ratio 

hrudndl=  hrudl/L  ;%  Non  Dimensionalized  Height 

Cdrudl  =  0.008  ;%  (input)  Rudder’s  Coef  of  Drag  (0.008) 

mrudl  =  (2*pi)/(l+2/arudl);%  (input)  Rudder's  Lift  Curve  Slope  Calc 

ProfAmd2  =  1 17.4  ;%  (input)  Rudder  Profile  Area  (1 17.4) 

PAnond2  =  ProfArud2/(L*T)  ;%  Non  Dimensionalized  Profile  Area 

hrud2  =  1 1.95  ;%  (input)  Rudder  Height 

arud2  =  hrud2A2/ProfArud2  ;%  (input)  Rud  Eff  Aspect  Ratio 

hrudnd2=  hrud2/L  ;%  Non  Dimensionalized  Height 

Cdrud2  =  0.008  ;%  (input)  Rudder's  Coefficient  of  Drag  (0.008) 

mrud2  =  (2*pi)/(l+2/arud2);%  (input)  Rudder’s  Lift  Curve  Slope  Calc 

ProfAdom  =  1400  ;%  (input)  Sonar  Dome  Profile  Area  (1400) 

PAnond  =  ProfAdom/(L*T)  ;%  Non  Dimensionalized  Profile  Area 

adorn  =  0.35  ;%  (input)  Sonar  Dome  Aspect  Ratio 

hdom  =  22.26  ;%  (input)  Sonar  Dome  Height 


47 


hdomnd=  hdom/L 


;%  Non  Dimensionalized  Height 


%  CALL  TO  INERTIAL  COEFFICIENT  GENERATOR 
csgen 


%  DERIVATIVE  GENERATOR  FOLLOWS  EQN'S  FROM  REF.  1  CHAP-8  pp525 
%  EQN’S  67  FOR  HULL  WITH  SMALL  DEADWOOD 

%  NON-DIMENSIONAL  HULL  HYDRODYNAMIC  DERIVATIVE  GENERATOR 
mpri  =  mship/(0.5*LA3*rho); 

Yvh  =  -pi*(T/L)*Cdo; 

Yrh  =  -(kl  *mpri)  +  ((xp/L)*Yvh); 

Nvh  =-(m2-(kl*mpri))  +  ((xp/L)*Yvh); 

Nrh  =  -(mz*(xbar/L))  +  ((.5*Priz)A2*Yvh); 

Yvdoth  =  -m2; 

Yrdoth  =  0; 

Nvdoth  =  0; 

Nrdoth  =  -((kprime*pi)/(LA3*TA2))*int3; 

Izh  =  mship/(0.5*LA5*rho); 

Kgpri  =  Kg*(L/U); 

Xgpri  =  Xg/L; 
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%  NON  DIMENSIONAL  RUDDER  (CONTROL)  HYDRODYNAMIC  DERIVATIVE 
GENERATOR 

Ydel  =(4/3)*PAnondl  *((2*pi)/(l+(2/arudl)));  %  4/3  ACCOUNTS  FOR  TWO  RUDDERS 
Ndel  =  -.5*(Ydel); 


%  NON  DIMENSIONAL  APPENDAGE  HYDRODYNAMIC  DERIVATIVE 

Generator 

Yvfirl  =  -PAnond  1  *((2*pi)/(  1  +(2/arud  1 ))); 

Yvfr2  =-PAnond2*((2*pi)/(l+(2/arud2))); 

Yvfd  =  -PAnond*((2*pi)/(l+(2/adom))); 

Yvap  =  Yvfrl+Yvff2+Yvfd; 

Yrfrl  =  -.5*Yvfrl; 

Yrff2  =  -.5*Yvfr2; 

Yrfd  =  -.5*Yvfd; 

Yrap  =  Yrfr  1  +  Y  rfr2+Y rfd; 

Nvfrl  =  Yrfrl; 

Nvfr2  =  Yrfr2; 

Nvfd  =  Yrfd; 

Nvap  =  Nvff  1+Nvff2+Nvfd; 

Nrfrl  =  .25*Yvfrl; 

Nrfr2  =  ,25*Yvfr2; 
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Nrfd  =.25*Yvfd; 


Nrap  =  Nrff  1+Nrfr2+Nrfd; 

Yvdotfrl  =  -4*pi*PAnondl*hrudndl/sqrt(arudl+l); 
Yvdotfr2  =  -4*pi*PAnond2*hrudnd2/sqrt(arud2+l); 
Yvdotfdm  =  -4*pi*PAnond*hdonmd/sqrt(adom+l); 
Yvdotap  =  Yvdotfrl+Yvdotfr2+Yvdotfdm; 

Yrdotfrl  =  2*pi*PAnondl*hrudndl/sqrt(arudl+l); 
Yrdotfir2  =  2*pi*PAnond2*hrudnd2/sqrt(arud2+l); 
Yrdotfdm  =  2*pi*PAnond*hdomnd/sqrt(adom+l); 
Yrdotap  =  Yrdotfrl+Yrdotfr2+Yrdotfdm; 

Nvdotfrl  =  Yrdotfrl; 

Nvdotfr2  =  Yrdotff2; 

Nvdotfdm  =  Yrdotfdm; 

Nvdotap  =  Nvdotff  l+Nvdotfr2+Nvdotfdm; 

Nrdotfrl  =  .25*  Yvdotfrl; 

Nrdotff2  =  ,25*Yvdotff2; 

Nrdotfdm  =  .25*Yvdotfdm; 

Nrdotap  =  Nrdotfrl +Nrdotfr2+Nrdotfdm; 

%  TOTAL  HYDRODYNAMIC  GENERATOR 
Yvtot  =  Yvh  +  Yvap; 
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Yrtot  =  Yrh  +  Yrap; 

Nvtot  =  Nvh  +  Nvap; 

Nrtot  =  Nrh  +  Nrap; 

Yvdot  =  Yvdoth  +  Yvdotap; 
Yrdot  =  Yrdoth  +  Yrdotap; 
Nvdot  =  Nvdoth  +  Nvdotap; 
Nrdot  =  Nrdoth  +Nrdotap; 

%  MATRIX  GENERATION 
al  1  =  mpri-Yvdot; 
al2  =  (mpri*Xgpri)- Yrdot; 
a21  =  (mpri*Xgpri)-Nvdot; 
a22  =  Izh-Nrdot; 

bll  =  Yvtot; 
bl2  =  Yrtot-mpri; 
b21  =  Nvtot; 

b22  =  Nrtot-(mpri*Xgpri); 

ell  =  Ydel; 
c21  =  Ndel; 
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A  =  [all  al2;a21  a22]; 
B  =  [bllbl2;b21b22]; 
C  =  [ell  ;c2 1  ] ; 


Apri  =  inv(A)*B; 

Bpri=(inv(A)*C); 

Cpri  =  [1  0;0-l]; 

Dpri=[0;0]; 

%  WIND  COEFFICIENT  GENERATOR  FOR  POOR  THESIS 
%  CONSTANTS 

rhow  =  64  ;%  Density  of  H20  in  lbs/ftA3 

rhoa  =  0.0752  ;%  Density  of  Air  in  lbs/ftA3 

LBP  =  390.00 ;%  (input)  Ship  Length  Between  Perpendiculars 

LOA  =  409.31 ;%  (input)  Ship  Length  Overall 

Af  =  2037.9  ;%  (input)  Ship  Frontal  Area 

As  =  10* Af  ;%  (input)  Ship  Side  Area  (Estimate) 

q=rhoa/rhow; 

%  COEFFICIENT  GENERATION 
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Cwx  =  q*Af/LBPA2*7000  ;%  Wind  Surge  Force  Coefficient 

Cwy  =  q*As/LBPA2*8000  ;%  Wind  Sway  Force  Coefficient 

Cwn  =  q*As*LOA/LBPA3*1000 ;%  Wind  Yaw  Moment  Coefficient 


%  ENGINE  CONSTANTS  AND  TORQUE  GENERATION 
1=1 .9e5  ;%  Prop  shaft  moment  of  inertia 

gm=(32.2A2/(5721  *2240))  ;%  Equivalent  form  of  g/m 

Qt=219597  ;%  Input  torque  at  16  knots 

Kl=2*gm*0.185*1.9914*15.5A4  ;%  Tp,  Prop  Thrust  Calc. 

K2=gm*171.7658  ;%  Rs,  Ship  Resistance  Calc. 

K3=(6000/(2*pi*I))  ;%  Qf,  Fric.  Torque  Loss  Calc. 

K4=(0.0417*1.9914*15.5A5/0.985)*(l/(2*pi*I)) ;%  Qp,  Prop.  Torque  Loss  Calc. 
K5=l/(2*pi*I)  ;%  Qt  Multiplier 

Aeng=[-2*K2*27.024  2*K1*1.72;0  (-2*K4*1.72)-K3J; 

Beng=[0  K5]'; 

Ceng=[l  0]; 

Deng=[0]; 

[numeng,deneng]=ss2tf(Aeng,Beng,Ceng,Deng); 
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%%%  torquefit 

Torq  =  [219597  586322  689169];%  Shaft  Torques  at  Vel,  Speeds  From 
Vel  =  [16  25.26  26.49]  ;%  RDS-2010  Final  Report 

Speeds  =  linspace( 16,22,60)  ;% 

Torques  =  interpl(Vel,Torq,Speeds,'linear');% 

%plot(Speeds,T  orques) 

vin  =  [16  17  18  19  20  21]; 

tout  =  interpl(Vel,Torq,vin,,linear'); 
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APPENDIX  B  CSGEN.M 


%  FILE  TO  GENERATE  kl,k2,kprime,  AND  Cs  FOR  RDS  2010 
kl=0.3*(2*T/L);  %  Empirical  Formulae 

k2=l-(.5*(2*T/L));  %  From  Vann  Thesis 

kprime=  1-(1.33*(2*T/L));  % 


%  SECTION  AREA  AT  EACH  STATION 

secar=[0  0  0  25.69  52.5  158.7  271.48  382.65  485.27  574.11... 

646.09  700.07  736.36  755.95  759.89  748.75  722.26... 

668.83  585.81  470.15  324.39  164.23  30.00];  . 

%  LOCAL  DRAFT  AT  EACH  STATION 

d=  [0  0  0  9.61 1  15.006  15.006  15.006  15.006  15.006  15.006... 

15.006  15.006  15.006  15.006  15.006  15.006  15.006  14.743... 

13.861  12.187  9.516  5.716  1.308]; 

%  LOCAL  BEAM  AT  EACH  STATION  (PAD  W/ 1  'S  TO  AVOID  DIVIDE  BY  0  ERROR) 
b=  [0  0  0  2.452  4.82  14.418  24.208  33.062  40.386  45.984  49.934... 

52.488  53.978  54.734  55.032  55.050  54.844  54.196  52.798  50.206... 

45.996  39.996  32.574]; 
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%  DISTANCE  x  FROM  CL  AT  EACH  STATION 
x=  [204.655  194.997  185.34  179.93  174.521  154.704  134.887  115.069... 
95.252  75.435  55.618  35.801  15.984  -3.834  -23.651  -43.468... 

-63.285  -86.848  -110.410  -130.973  -157.535  -181.098  -204.66]; 

for  i=l:length(d) 
if  d(i)  =  0 
mult(i)  =  0; 
else 

mult(i)  =  secar(i)/(b(i)  *d(i)) ; 
end 

if  b(i)  =  0 
mult2(i)  =0; 
else 

mult2(i)  =  (4*d(i)/b(i)); 
end 

co(i)  =  -0.8572  +  0.5339*mult2(i);  %  FROM  FORMULAE  IN  REF.  3 

cl(i)  =  3.734  -  1.3661*mult2(i); 
c2(i)  =  -1.7323  +  0.8679*mult2(i); 

Cs(i)  =  co(i)+(cl(i)*mult(i))+(c2(i)*mult(i)A2); 
if  Cs(i)  >  0.918 
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Cs(i)  =  0.918 


;%  CAPS  MAX  Cs  VALUE  AT  MAX  AREA  STA. 


end 

Csdsq(i)  =  Cs(i)*  d(i)A2;  %  CALCULATES  VECTORS  FOR  INTEGRATION 

Csdsqx(i)  =Cs(i)*d(i)A2*x(i); 

Csdsqxsq(i)=Cs(i)  *d(i)A2*x(i)A2 ; 
end 

%  DISTANCE  OF  EACH  ALONG  HULL  FROM  BOW  TO  STERN 
huldist=  [0  9.658  19.315  24.725  30.134  49.91  69.768  89.586  109.403... 

129.22  149.037  168.854  188.671  208.489  228.306  248.123... 

267.94  291.502  315.065  335.628  362.19  385.753  409.315]; 
xi=linspace(0,409,409) ; 


Csdsqfit=interp  1  (huldist,Csdsq,xi)  ;%  FIT  STATION  DATA  TO  LENGTH  OF  HULL 

Csdsqxfit=interp  1  (huldist,Csdsqx,xi)  ;%  FIT  STATION  DATA  TO  LENGTH  OF  HULL 

Csdsqxsqfit=interpl(huldist,Csdsqxsq,xi);%  FIT  STATION  DATA  TO  LENGTH  OF  HULL 

int  1  =trapz(Csdsqfit) ; 

int2=trapz(Csdsqxfit) ; 

int3=trapz(Csdsqxsqfit); 

m2=  (k2*pi)/(L*TA2)*intl; 

xbar=(int2)/intl; 


57 


mz=(kprime/k2)*m2; 


kl; 

k2; 

kprime; 

%subplot(21  l),plot(Cs) 
%subplot(212),plot(Csfit) 
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